module paint_on_mod

use global_mod

use domain_mod,     only: Var_Find, Var_Get
use mesh_mod,       only: MeshObj
use utils_misc_mod, only: Critical_Error, ARROW_ERR_PON_FORM, ARROW_ERR_PON_FUNC

implicit none

character(len=32) :: form     = "cube"
character(len=32) :: function = "flat"
character(len=16) :: variable = ""
real(kind=rlk)    :: A        = 1.0
real(kind=rlk)    :: x0(3)    = (/0.0,0.0,0.0/) !< form centre
real(kind=rlk)    :: dx(3)    = (/1.0,1.0,1.0/) !< form width
real(kind=rlk)    :: f0(3)    = (/0.0,0.0,0.0/) !< function centre
real(kind=rlk)    :: df(3)    = (/1.0,1.0,1.0/) !< function width
namelist /painton/ form, function, variable, A, x0, dx, f0, df

integer(kind=ink), parameter :: ARROW_FORM_CUBE   = 1
integer(kind=ink), parameter :: ARROW_FORM_SPHERE = 2
integer(kind=ink), parameter :: ARROW_FORM_CONE   = 3
integer(kind=ink), parameter :: ARROW_FUNC_FLAT   = 1
integer(kind=ink), parameter :: ARROW_FUNC_GAUSS  = 2

type PaintObj
  integer(kind=ink) :: form
  integer(kind=ink) :: func
  integer(kind=ink) :: ivar 
  real(kind=rlk)    :: A
  real(kind=rlk)    :: x0(3)
  real(kind=rlk)    :: dx(3)
  real(kind=rlk)    :: f0(3)
  real(kind=rlk)    :: df(3)
end type

type(PaintObj), pointer :: PaintList(:)
integer(kind=ink)       :: nPaint

contains


type(PaintObj) function Paint_Create() result(pon)

  select case(form)
    case("cube")
      pon%form = ARROW_FORM_CUBE
    !case("sphere")
    !  pon%form = ARROW_FORM_SPHERE
    !case("cone")
    !  pon%form = ARROW_FORM_CONE
    case default
      call Critical_Error(ARROW_ERR_PON_FORM,"PaintCreate")
  end select

  select case(function)
    !case("flat")
    !  pon%func = ARROW_FUNC_FLAT
    case("gaussian")
      pon%func = ARROW_FUNC_GAUSS
    case default
      call Critical_Error(ARROW_ERR_PON_FUNC,"PaintCreate")
  end select

  pon%ivar = Var_Find(variable)
  
  pon%A  = A
  pon%x0 = x0
  pon%dx = dx
  pon%f0 = f0
  pon%df = df
end function Paint_Create


subroutine Paint_Setup(nn)
  integer(kind=ink), intent(in) :: nn
  allocate(PaintList(nn))
  nPaint = nn
end subroutine Paint_Setup


subroutine Paint_On(mesh)
  type(MeshObj), intent(inout) :: mesh
  type(PaintObj), pointer :: pon
  integer(kind=ink) :: ii, pp
  real(kind=rlk), pointer :: var(:)
  real(kind=rlk) :: xx(1:3)
  do pp=1,nPaint
    pon => PaintList(pp)
    var => Var_Get(pon%ivar)
    do ii=1,mesh%nCells
      !xx = mesh%cell_x(ii)
    enddo
  enddo
end subroutine Paint_On


real(kind=rlk) function Cube(xx, x0, dx) result(ff)
  real(kind=rlk), intent(in) :: xx(3)
  real(kind=rlk), intent(in) :: x0(3)
  real(kind=rlk), intent(in) :: dx(3)
  integer(kind=ink) :: ii
  ff = 1.0
  do ii=1,3
    if (xx(ii)<x0(ii)+dx(ii)) ff = 0.0
    if (xx(ii)>x0(ii)-dx(ii)) ff = 0.0
  enddo
end function Cube


real(kind=rlk) function Gaussian(xx, x0, dx) result(ff)
  real(kind=rlk), intent(in) :: xx(3)
  real(kind=rlk), intent(in) :: x0(3)
  real(kind=rlk), intent(in) :: dx(3)
  integer(kind=ink) :: ii
  ff = 1.0
  do ii=1,3
    ff = ff*exp(-(xx(ii)-x0(ii))**2/dx(ii)**2)
  enddo
end function Gaussian

end module paint_on_mod